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The propagation of filaments in the Scrape Off Layer (SOL) of tokamaks largely determine the plasma profiles 
in the region. In a conduction limited SOL, parallel temperature gradients are expected, such that the 
resistance to parallel currents is greater at the target than further upstream. Since the perpendicular motion 
of an isolated filament is largely determined by balance of currents that flow through it, this may be expected 
to affect filament transport. 3D simulations have thus been used to study the influence of enhanced parallel 
resistivity on the dynamics of filaments. Filaments with the smallest perpendicular length scales, which were 
inertially limited at low resistivity (meaning that polarization rather than parallel currents determine their 
radial velocities), were unaffected by resistivity. For larger filaments, faster velocities were produced at higher 
resistivities, due to two mechanisms. Firstly parallel currents were reduced and polarization currents were 
enhanced, meaning that the inertial regime extended to larger filaments, and secondly a potential difference 
formed along the parallel direction so that higher potentials were produced in the region of the filament for 
the same amount of current to flow into the sheath. These results indicate that broader SOL profiles could 
be produced at higher resistivities. 


I. INTRODUCTION 

One of the biggest challenges facing future generations 
of magnetic confinement devices such as ITER is to con¬ 
trol the high particle and heat fluxes at the divertor sur¬ 
faces. These fluxes are determined by the balance be¬ 
tween transport across and parallel to the magnetic field 
in the Scrape Off Layer (SOL), with enhanced cross field 
transport leading to a broader SOL width and hence re¬ 
duced fluxes to the divertor. 

Perpendicular transport of particles (and to a lesser 
extent heati) in the SOL has been observed to be 
dominated by the radial advection of coherent plasma 
structures^, that are significantly more dense and hot 
than their surrounding plasma, with peak fluctuations 
typically of the order of the background^. These struc¬ 
tures are aligned to the equilibrium magnetic field and 
strongly localized in the drift-plane perpendicular to it 
and hence are referred to as filaments (or hlohs due to 
their appearance in the perpendicular plane). A num¬ 
ber of recent works have provided reviews of the exper¬ 
imental evidence for filaments^, the theoretical under¬ 
standing of their propagation and contribution to SOL 
transport^i^, and of the agreement found when exper¬ 
imental measurements of filaments are compared with 
theory and simulations^. 

The basic mechanism by which filaments advect ra¬ 
dially outwards from the core was first proposed by 
Krasheninnikovi and can be understood in a fluid model^ 
by considering the balance of electrical currents through 
an isolated filament in the SOL. In the region of the out- 
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board mid-plane of a tokamak, the gradient and curva¬ 
ture of the magnetic field act to enhance diamagnetic 
currents at larger radial distances for a given pressure 
gradient. The cross-field pressure gradients within the 
filament thus lead to diamagnetic currents in the per¬ 
pendicular plane that have a non-zero divergence. Cur¬ 
rent continuity necessitates additional currents, which 
can take the form of perpendicular polarization currents 
to produce a closed circuit within the drift-plane, or par¬ 
allel currents that can close through the sheath at the 
target. Each of these current paths result in the forma¬ 
tion of a broadly dipolar electrostatic potential field in 
the perpendicular plane, which through E x B motions 
corresponds to a pair of counter rotating vortices that 
act to advect the filament radially outwards. The mag¬ 
nitude and structure of this potential field and thus the 
filament’s detailed motion, is dependent on whether the 
parallel or polarization current path is dominant in clos¬ 
ing the diamagnetic currents, which is in turn determined 
by the effective resistances of each path. 

Much of the work to date concerning the theory and 
simulation of filament motion has concentrated on the 
case in which the resistance to currents traveling through 
the sheath, Tsheath, is much greater than the total resis¬ 
tance to parallel currents through the bulk SOL plasma, 
Fy, so that resistance of the entire parallel current path 
is dominated by T sheath- In this limit, it has been shown 
that the perpendicular length scale of the filament per¬ 
turbation, i5_l, plays an important role in determining 
whether the parallel or polarization current path is dom¬ 
inant in closing the diamagnetic current drive, and thus 
the filament’s motions^d^. Estimates of how a filament’s 
radial velocity scales with this parameter have been de¬ 
rived. For filaments much smaller than a critical length, 
i5_l <C <5*0, where <5*o is defined in Equation m, the po- 
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larization current path is dominant and the filament’s ve¬ 
locity is estimated to scale like ~ Such filaments 

are said to be in the inertial regime or inertially lim¬ 
ited. On the other hand, for filaments much larger than 
this critical length, i5_l ^ (5*o, parallel currents traveling 
through the sheath to close at the target become domi¬ 
nant and the velocity is estimated to scale like ^ In 

this case, the filaments are described to be in the sheath 
current regime or sheath current limited. These asymp¬ 
totic regimes have been observed in two dimensional^!^, 
and more recently three dimensional, simulation s^^ d ^ . 

There exist a number of mechanisms however by which 
r |j may significantly increase, such that the total resis¬ 
tance of the entire parallel current path is no longer 
dominated by Tsheath- One such mechanism is an in¬ 
creased parallel connection length to the target, which 
can be achieved in a larger device or by use of a Super-X 
divertor—. Moreover, decreased temperatures within the 
SOL will lead to an enhanced electron-ion collisionality 
and thus a higher parallel resistivity. In a conduction lim¬ 
ited SOL this effect can become very strong particularly 
in the divertor region, as the temperature downstream 
at the target can become much cooler than upstream at 
the mid-plane. Moreover, if the temperature at the tar¬ 
get becomes sufficiently low (Te <leV), volume recom¬ 
bination becomes strong and divertor detachment will 
occur—, meaning that a cloud of neutrals forms between 
the plasma and the target. Once the ionization fraction 
is sufficiently small, electron-neutral collisions are com¬ 
parable with electron-ion collisionsi^, and the resistivity 
has a component proportional to the ratio between neu¬ 
tral and electron densitiesi^ (see appendix). In the limit 
of zero ionization, the resistance of the neutral gas in 
front of the targets is effectively infinite. 

The effect of parallel resistivity was considered using a 
two region model in Reference [l3- In the absence of mag¬ 
netic geometry effects which were also considered, the 
work’s predictions are equivalent to the aforementioned 
inertial and sheath current regimes when Tsheath ^ ry 
(although they are referred to as the resistive ballooning 
and sheath connected interchange regimes respectively). 
In the case in which ry > Tsheath, the inertial regime is 
expected to continue for the smallest <5_l, as the regime 
does not involve parallel currents. For larger dj_ filaments 
in this case however, what is described in this work as a 
resistive sheath current regime is predicted (Reference 
[TtI uses the term resistive X-point regime), in which the 
diamagnetic currents are closed through the sheath, but 
the parallel resistance of the plasma, rather than that of 
the sheath is expected to dominate in determining the 
electrostatic potential within the filament, and thus its 
radial velocity. The radial velocity within this regime is 
predicted to scale like ~ r||(5j_“^. This two region model 
also predicts that the critical d_L at which filaments tran¬ 
sition from the inertial to the resistive sheath current 
regime scales like <5* ~ Therefore if the plasma’s 

parallel resistivity is sufficiently high, then effectively all 
filaments of a physically realistic d_L will be in an inertial 


regime. 

It is this extreme limiting case that the two dimen¬ 
sional ESEL model— considers by neglecting the influ¬ 
ence of parallel currents entirely. Turbulence simulations 
using this model have been successful in reproducing ex¬ 
perimental measurements of SOL profiles and turbulence 
statistics from a variety of experimental devicesi^^— . 
The absence of parallel current effects within this model 
mean that it only reproduces the inertial regime. Iso¬ 
lated filament simulations using the ESEL model there¬ 
fore did not find good a gree ment with three dimensional 
simulations in Reference I 111, which used parameters such 
that Tsheath ^ FII and so also exhibited the sheath cur¬ 
rent regime. It was however suggested that inclusion 
of additional physics such as an enhanced collisionali- 
ties or divertor detachment within the three dimensional 
model would make its results more comparable to those 
obtained using the ESEL model. 

The three dimensional simulations presented in this 
work therefore investigate the effect of enhanced resis¬ 
tivity, particularly in the region nearest to the targets, 
on the current balance within filaments and hence their 
radial velocity. The remainder of this paper is organized 
as follows. Section [IT] provides an outline of the physical 
model used in this work, before Section |nT] describes the 
numerics, boundary conditions and initialization of the 
simulations. The results in Section HVl are split into three 
parts. Eirst Section HV Al describes the dynamics of a low 
resistivity case that is used as a reference case which the 
enhanced resistivity simulations are compared against. 
Next, the effect of enhancing the resistivity only in the 
region nearest the target is investigated in Section lIVBi 
before the Section IIV Cl demonstrates the effect of resis¬ 
tivity enhanced uniformly throughout the entire domain. 
Einally the main conclusions of this work are summarized 
in Section |Vl 


II. PHYSICAL MODEL 

The simulations presented in this paper have been ob¬ 
tained using the same physical model as in Reference fill, 
which is an electrostatic drift-fluid model that assumes 
singly charged cold ions and isothermal electrons. It is 
acknowledged that the assumption of cold ions is poorly 
justified in the SOL, where typically Ti > but it is 
used in this work for simplicity. 

The effects of magnetic geometry have been neglected 
by employing a slab geometry with uniform magnetic 
field B = Bz to represent the SOL, with the effects of 
magnetic curvature and gradients included through ad¬ 
ditional terms in the evolution equations. The radial 
coordinate in this geometry is represented by the x coor¬ 
dinate, whilst y corresponds to an effective poloidal co¬ 
ordinate. Throughout this work a Bohm normalization 
is used, with time and length scales normalized to the 
ion gyro-frequency, 17^ = eR/rUi, and the hybrid gyro- 
radius ps = Cs/^i respectively, whilst the electrostatic 








3 


potential, (j), has been normalized to Tg/e. Here e is the 
elementary unit charge, rui is the ion mass, Cs = ^/Te/mg 
is the sound speed, Tg is the electron temperature in 
Joules and rUg is the mass of an electron. In addition, 
the plasma density has been normalized to a character¬ 
istic SOL plasma density, ng. The resulting dimension¬ 
less evolution equations for plasma density, n, vorticity, 
O = parallel ion velocity, t/, and parallel electron 

velocity, V, are: 


f = + + (1) 

_ =-V|| (nI/)-hng^-g—S'™, (2) 


dt 


H' 


SnU 

n 


(3) 


^ J|| - (4) 


Here, ^ ^ x ' V), J|| = n{U — V) is the nor¬ 

malized parallel current density, S„ is a particle source, 
/i = mil mg is the ratio of ion to electron masses, Dn is 
the normalized particle perpendicular diffusivity, /i^ is the 
normalized ion perpendicular viscosity, = z/ezo/l-96Hi, 

I'eio = noe^ In A/3my^eQ(27rre)^/^ is the electron-ion col¬ 
lision frequency, In A is the Coulomb logarithm and eg is 
the permittivity of free space. It is through the param¬ 
eter j^ii that the effective parallel resistivity, ? 7 || = 
was controlled for the studies in this paper. 

Equation 0 enforces current continuity and is simply 
the divergence of current density divided through by the 
plasma density, V ■ J /n, under the Boussinesq approxi¬ 
mation. Whilst the second term on the Right Hand Side 
(RHS) of this equation is written explicitly as V|| Jy/n, 
the remaining terms within this equation originate from 
the other currents in the system as follows. The left hand 
side and first term on the RHS result from the ion polar¬ 
ization current, Jpoi- The electron diamagnetic current 
density, Jdia = bx VnjB leads to the third term on the 
RHS, with the last term arising from the viscous current 
density, Jvisc, which exists due to the presence of viscos¬ 
ity in the system. In other works, Jvisc is often included 
within Jpoi- It is through the diamagnetic current that 
the VH and curvature effects force the other currents in 
the system. The strength of such a drive has been repre¬ 
sented through the dimensionless variable g, which at the 
outboard mid-plane of a tokamak can be approximated 
to be g = 2psjRc, where Rg is the dimensional radius 
of curvature. Concerning the other governing equations 
of the model. Equations © to dH) conserve particle den¬ 
sity and parallel momentum for each particle species with 
Equations ([3]) and ([3]) written in non conservation form. 

At the location of the entrance to the sheath in front of 
the target, z = ±L||, where Ly = i\\/ps is the normalized 
mid-plane to target distance, the parallel velocity fields 
(which are normalized Cg) evolved by these two equations 


must satisfy standard sheath boundary conditions'^: 



> 1, 

(5) 


< -1, 

(6) 

-H 

II 

= ±exp(-((>|^^±^J . 

(7) 


In writing equation © in such form, (j) is dehned relative 
to the potential at the target plate wall <j)yj , which is given 
a fixed value: 


(j)ni = - In 



( 8 ) 


The target is thus assumed to have zero resistance. Erom 
Equations ([5]) to 0, the resistance to parallel currents 
traveling through the sheath and closing at the target 
can be estimated. Assuming the potential at the sheath 
entrance to be small, Equation © can be linearized so 
that the current density traveling through the sheath is 
approximately Jy g ~ ns(j)s, where the subscript s de¬ 
notes the value of a field at the entrance to the sheath. 
Two points at the sheath entrance which have a potential 
difference with equal magnitude but opposite sign corre¬ 
spond therefore to the parallel current density at the first 
point traveling through the sheath, closing through the 
target and emerging at the second point. The potential 
difference between the two points is 2J||/ng and by us¬ 
ing basic circuit theory the resistance (for a unit area) 
to currents traveling through the sheath to close at the 
target is 2/ng. Since this path involves the current going 
through the sheath twice, the effective resistance of the 
sheath itself is: 

^sheath — ■ (9) 

Us 

The effective resistance of parallel currents traveling from 
the diamagnetic current source, assumed to be located at 
the mid-plane {z = 0), to the sheath entrance can also 
be calculated: 



( 10 ) 


Using these effective resistance definitions, the theoret¬ 
ical estimates from References and 113 for the critical 
value of i5_l at which filaments transition from the inertial 
regime, S^,, can be summarized as: 


where 


J d*g Ty <C Tg heath 

\ d*gry^/® Ty » Tg/jeat^i 



( 11 ) 


( 12 ) 


The Ty <C Tsheath condition corresponds to the transition 
between the inertial and sheath current regimes, whilst 
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the rII ^ Tsheath condition corresponds to Reference H^’s 
prediction for a transition from the inertial regime to 
the anticipated resistive sheath current regime. Equation 
(1121) can be derived by following the scaling arguments in 
References d or [2^ It is noted from the former reference 
that (5*0 should also include an order unity correction to 
account for the magnitude of the density perturbation, 
but since this can only be determined numerically via an 
amplitude scan, it has been neglected here. 


III. SIMULATION IMPLEMENTATION 

A. Numerics 

The results within this paper were obtained us¬ 
ing the SOL2Fluid physics moduleii written using the 
BOUT-I--I- frameworb22i^. The time integration was car¬ 
ried out using a fully implicit Newton-Krylov Backwards 
Difference Formula (BDF) solver from the PVODE li¬ 
brary. All the spatial derivatives were calculated using 
second order accurate schemes and thus the code was sec¬ 
ond order accurate. Specifically, the parallel advection 
derivatives were calculated using an upwind scheme, an 
Arakawa scheme^^ was used for the perpendicular E x B 
advective terms and all other derivatives were calculated 
using central differencing. For numerical stability, the U 
and V fields were staggered in the parallel direction rel¬ 
ative to the other fields. Both the BOUT-I--I- framework 
and the SOL2Fluid physics module has been successfully 
verified using the Method of Manufactured Solutions^^. 



FIG. 1. Variation along the parallel direction of the equilib¬ 
rium fields used for all of the simulations in this work. The 
fields are all uniform in the perpendicular plane. 


C. Initialization 


B. Boundary Conditions 

For computational efficiency, only half the parallel do¬ 
main was simulated, with symmetry boundary conditions 
employed at the lower parallel boundary at z = 0. At 
the upper parallel boundary, z = L||, Equations (O and 
© were enforced on U and V. No boundary conditions 
were specified for the remaining variables n, ipj ^ the 
upper parallel boundary to avoid over constraint of the 
system. In the perpendicular plane, the y direction was 
periodic for all fields, whilst at the x boundaries Q and 
the gradients of n, U and V were set to zero. The re¬ 
maining perpendicular boundary condition for (j) was set 
to obtain the ID (variation only in the parallel direction) 
equilibrium fields described in the next subsection. This 
was achieved by fixing the x boundaries of (j) to the par¬ 
allel profile of its resulting equilibrium field. However, as 
this boundary condition could not be determined a pri¬ 
ori^ it was obtained by iteratively running the simulation 
until it achieved a steady state equilibrium and updating 
the (j) boundary condition to its parallel profile along the 
center of the domain until a time invariant ID system 
was produced. 


In order to isolate the dynamics of single filaments from 
that of the background fields, a steady state equilibrium 
with variation only in the parallel direction was required, 
onto which the filament density perturbations could be 
seeded. This was achieved by evolving the system until 
time invariant background fields were obtained using the 
following density source: 


10 exp (l0z/L||) 
L\\ (exp (10) - !)■ 


(13) 


The source is predominantly localized in the last 10% 
of the domain nearest the target, and this structure was 
chosen as it produces equilibrium fields with negligible 
parallel velocities and parallel gradients of density and 
potential for the majority of the domain. The source 
structure can also be interpreted to loosely model a high 
recycling regime. The equilibrium field profiles (denoted 
by the subscript eq) produced using this source are shown 
in Figure [U with the exception of ileq, which is necessar¬ 
ily zero. This equilibrium has been verified against ana¬ 
lytical resultsii to ensure that the equations were being 
solved correctly. 

The filaments were modeled within the simulations as 
density perturbations on top of the equilibrium density 
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field. It is emphasized that the simulation code evolved and that the equilibrium and filament perturbation were 
the quantities according to the full non linear equations, not evolved separately. Each simulation was seeded with 

a single density perturbation of the form 

I 


I O I ( x'^+y^ 


where rif = n—neq, is the density perturbation of the fila¬ 
ment. The filaments were thus seeded as Gaussian struc¬ 
tures in the perpendicular plane with a perpendicular 
length scale dx- In the parallel direction they extended 
approximately half the length of the domain from the 
mid-plane at z = 0 to z = L||/2. This structure along the 
field line was chosen as filaments are believed to be gen¬ 
erated inside of the separatrix in the outboard mid-plane 
regional, and have been observed in MAST to extend 
from X-point to X-point within the SOL^. Moreover, 
3D filament simulations that included realistic magnetic 
geometry have shown that filaments initialized extending 
from target to target will rapidly develop parallel density 
gradients due to ballooning motions around the outboard 
mid-plane and enhanced dissipative effects downstream 
of the X-point region^Si. 

Each simulation utilized a grid mesh of x Ny xNz = 
192 X 128 X 16 grid points, where Ni is the number of grid 
points in the i coordinate direction, and the perpendic¬ 
ular domain size was scaled with (5 _l so that the lengths 
of the simulation domain in the x and y directions were 
Lx = 15dj_ and Ly = 10dj_ respectively. 
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target at z = L||, the density at the sheath remained 
largely constant at its equilibrium value, Ug ss 0.85 (See 
Figure [T]) and so Tsheath ~ 1-2 all of the simulations 
presented in this work. 


A. Reference Case 

In order to demonstrate the effect of resistivity, the dy¬ 
namics of filaments using the parameters given in Table 
H] are hrst described. The results using these parame¬ 
ters, which are such that Tji -C Lsheath, will then be used 
as a reference case against which the simulations with 
higher resistivity will be compared. The parameters are 
broadly relevant to the conditions found in the Mega Am¬ 
pere Spherical Tokamak (MAST^ and were chosen to 
allow comparison with Reference I 111 The dissipative pa¬ 
rameters Dn and yt were specified to be two orders of 
magnitude smaller than neoclassical estimates for their 
valuesi^ (and so in effect broadly classical^! values were 
used) to ensure that viscous currents played a negligible 
role. 


IV. RESULTS 

Each simulation was evolved until t = 500 with an out¬ 
put time step of 5 for dj_ < 5 and 25 for (5j_ > 5, which 
was sufficient to obtain a characteristic radial velocity, 
u/, of the filament. This characteristic velocity was de¬ 
fined as the first maximum that occurred of which 
is the instantaneous radial center of mass velocity of the 
filament perturbation: 

OO OO OO 

A A A 

n f ax ay az 
dy 

— OO —OO —OO 

OO OO OO ' 

nf dxdy dz 

— OO —OO —OO 

It is noted that was not monotonic, and the first 
maximum was selected to avoid the effects of Boltzmann 
spinning motion o^^'^^ that could produce a second larger 
maximum for the smallest (5j_ filaments. 

Since the simulations were not evolved for long enough 
for the density perturbation of the filament to reach the 


TABLE I. Reference Case Parameters 


Input Parameters 

Dimensionless Parameters 

Te = 40 eV 

g = 2.43 X lO""" 

B = 0.5 T 

nil = 2.53 X 10”^ 

no = 0.8 X lO^^cm'^ 

Dn = 1.8 X 10"® 

Rc — 1.5 m 

Pi = 5.0 X lO""* 

I'll = 10 m 

p = 3646 

mi = 3.32 X 10”^^ kg 

L|| = 5500 

InA = 13.3 

5*0 = 8.2 


Pll = 3.8 X 10"^ 


The dependence of u/ on for the reference set of 
parameters is shown in Figure [2] and it is clear that these 
simulations show good agreement with the T || <C Tsheath 
analytical scaling estimates discussed in Section HI which 
are also plotted for comparison. Filaments initialized 
with d_L <C (5*0 1 which is plotted using a grey dotted hor¬ 
izontal line, produced characteristic velocities that scale 
like (5j_^^^. On the other hand characteristic velocities 
proportional to i5j_~^ were obtained by hlaments initial¬ 
ized with (5 _l i5*o- For reference, the peak value of 

Vf = 0.46 corresponds to a dimensional radial velocity of 












0 


4 


6 


0.1 
0.09 
0.08 
0.07 
0.06 
^ 0.05 
0.04 
0.03 
0.02 
0.01 
0 

0 10 20 30 40 50 

FIG. 2. Dependence of the characteristic radial velocity, Vf, 
on its initial perpendicular length scale dx. The analytical 
scaling estimates for the inertial and sheath current regimes 
are also plotted for comparison. 
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2km/s, which is consistent with experimental measure¬ 
ments from MAST—. 

It is next demonstrated that the current balance found 
in each regime and thus the mechanisms by which the ve¬ 
locities are limited are also consistent with the theoretical 
predictions. The typical current balance found in large 
filaments, dx /§> <5*0; is displayed in Figure [31 which plots 
the divergences of each of the current densities divided 
by n (corresponding to the terms in Equation ([I} as de¬ 
scribed in Section m within the dx = 28 filament, in 
various drift-planes along the field line. The V • Jmscjn 
quantity has not been plotted as it is negligible. As in 
all subsequent contour plots of any quantity in this work, 
the quantities in this figure were taken at time at which 
the Vj occurred. The diamagnetic current drive can be 
seen to exist in the region from z = 0 to z = T||/2, and is 
almost entirely balanced by the parallel currents. That 
these parallel currents travel through the sheath and thus 
the filament can accurately be described as sheath cur¬ 
rent limited, is illustrated by the left plot of Figure [H 
which shows the structure of J|| in a ?/ — z plane through 
the middle of the filament, whose rif contours are over¬ 
laid for reference. The parallel currents, which are of 
opposite sign left and right of y = 0, can be seen to in¬ 
crease in magnitude along z from J|j = 0 at z = 0 until 
the end of the density perturbation at z ~ A||/2. They 
then remain relatively constant from there onwards as 
they travel through the background to the pre-sheath 
entrance at z = L||, where they can be interpreted to be 
traveling through the sheath to close through the target. 
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FIG. 3. Divergence of each current density divided by particle 
density, from the reference case (5x = 28 filament at various 
perpendicular planes along the field line. The quantities are 
plotted at the time at which the filament’s characteristic ra¬ 
dial velocity, Vf, occurred. 
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FIG. 4. Structure of Jy in a y — z plane through the middle of 
5x = 28 filaments at the time of their peak radial velocities. 
The left, center and right plots respectively correspond to 
the reference case, the enhanced target localized resistivity 
t^Wdiv = 10000iZ||o case, and the uniformly enhanced resistivity 
iz|| = 2500iZ||o case filaments. Each case is plotted at the 
time at which the filament’s characteristic radial velocity, Vf, 
occurred. 
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A more complicated current balance is typically ex¬ 
hibited in small filaments, <5^ <C ^*oj cis can be seen in 
Figure m which plots the divergences of current densities 
divided by n from the (5 _l = 5 filament. It is not imme¬ 
diately obvious from the top three rows of this Figure 
whether the polarization or parallel currents are play¬ 
ing a greater role in closing the diamagnetic currents. 
More clarity is provided by isolating the components of 
V • Jpoi/n and V • Jy/n which have even (-I-) and odd (—) 
parity in y with respect to the center of mass of the fila¬ 
ment in the poloidal direction, j/q- These decompositions 
can be calculated for an arbitrary quantity / as follows: 


where 


f f{y-yo)±f {yo - y) 

J± = -;;- 


oo oo oo 


J j J Ufy dx dy dz 


yo = 


— oo —oo —oo 
oo oo oo 


J J J Hf dxdydz 


— OO —OO —OO 


(16) 


(17) 


These decompositions are plotted in the bottom four rows 
of the figure. The diamagnetic currents can be seen to be 
predominantly balanced by the odd component of Jpoi 
and thus the filament’s velocity can be said to be iner- 
tially limited. An independent current balance is found 
between the even components of Jy and Jpou which arise 
as a result of the Boltzmann potential response to the 
parallel density gradient a^^i^^'^^ . The odd component of 
Jy also plays an important role however, as it balances 
the odd component of Jpoi in the region beyond z = Ty/2 
where the diamagnetic drive is no longer present. This 
means that the odd component of Jpoi and thus (j) is 
broadly constant all along the field line. 


B. Increased Target Localized Resistivity 


As discussed in SectionUl there are a number of mecha¬ 
nisms that may enhance the collisionality particularly in 
the divertor region. Investigations were therefore carried 
out to determine the effect of increasing the resistivity of 
the plasma in this region on the dynamics of filaments 
by increasing z/y from its reference case value, z/yo, in the 
last 25% of the domain nearest the target: 


f Z/yo Z<3Ly/4 

\ V\\div Z > 3Ly/4 


(18) 


Simulations were carried out using values of v\\div that 
increased from z/yo by powers of ten to lOOOOz^yo, cor¬ 
responding to values of Fy ranging from 0.038 to ap¬ 
proximately 95. For comparison Fy ss Tsheath when 
^\\div = An order of magnitude estimation of 

the electron temperature required in the divertor region 
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FIG. 5. Divergence of each current density divided by parti¬ 
cle density, from the reference case (5 _l = 5 filament at various 
perpendicular planes along the field line. The quantities are 
plotted at the time at which the filament’s characteristic ra¬ 
dial velocity, «/, occurred. 


to produce such enhancements in v^divi is provided in the 
appendix. 

The dependence of f/ on (5^ for increasing values of 
^\\div can be seen in Figure [S] and it is clear to see that 
increasing u^div leads to enhanced radial velocities across 
all 5j_, with the smallest experiencing a relatively 
modest increase in rzj, compared to the larger (5j_. This 
is to be expected because the smallest filaments were in 
the inertial regime in the reference case, meaning that 
parallel currents played a sub-dominant role in closing 
the diamagnetic currents. On the other hand, the largest 
hlaments were in the sheath current regime in the refer¬ 
ence case, meaning that parallel currents were dominant 
in maintaining current continuity and so increasing the 
resistivity has a greater influence on these filaments. 

One of the mechanisms by which faster velocities are 
produced can be understood by observing that as v\\div 
is increased, Vf scales like up until larger values 

of (5j_ and so the inertial regime is clearly extended, or 
equivalently is increased. This occurs because increas- 
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FIG. 6. Dependence of the characteristic radial velocity, Vf, 
on its initial perpendicular length scale dx for each of the 
values of v^idiv used in the enhanced target localized resistivity 
scan. 


ing the resistivity suppresses the parallel currents and 
thus necessarily leads to an enhancement of the polariza¬ 
tion currents, given the same diamagnetic current source. 
The drastic reduction of Jy in a (5x = 28 filament is 
evident upon comparison of the left and middle plots 
of Figure HI which plot this quantity in a, y — z plane 
through the center of the filament, for the reference case 
^\\div = i^llo and iy\\div = 10000z/||o simulations respec¬ 
tively. To demonstrate that this suppression of Jy affects 
the current balance upstream. Figure [7] plots the diver¬ 
gence of current densities divided by density, from the 
same simulation at the time of its peak radial velocity. 
By cross comparison with Figure [H it is clear that the 
reduction of parallel currents within the filament means 
that the diamagnetic current drive is predominantly bal¬ 
anced by enhanced polarization currents instead. The 
parallel currents do still play an important role however, 
in that they balance Jpoi in the region between the end of 
the filament perturbation and the start of the enhanced 
resistivity region, T||/2 < x < 3L||/4, where the diamag¬ 
netic drive is no longer present, so that Jpoi, and thus 
(j), is approximately constant in z up until the region of 
enhanced resistivity. This is analogous to the role they 
played for Jx ^ filaments in the reference case. 

Greater radial velocities were also produced at higher 
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FIG. 7. Divergence of each current density divided by par¬ 
ticle density, from the enhanced target localized resistivity 
5x = 28, = 10000iZ||o filament at various perpendicular 

planes along the field line. The quantities are plotted at the 
time at which the filament’s characteristic radial velocity, Vf, 
occurred. 


resistivities by the very largest i5x, in which parallel cur¬ 
rents closing at the target were still the main way in 
which the diamagnetic currents were closed. These fila¬ 
ments attained greater velocities because the resistance 
of the plasma was sufficient to introduce a potential dif¬ 
ference between the downstream at the sheath entrance 
and further upstream in the region of the filament den¬ 
sity perturbation. Therefore for the same amount of cur¬ 
rent to flow into or out of the sheath, larger potentials 
were formed upstream at higher resistivities, which in 
turn correspond to faster radial velocities. Such filaments 
are hereafter described to be in a resistive sheath current 
regime. 

To demonstrate the potential difference formed along 
the parallel direction at high resistivities, it is necessary 
to separate it from the potential difference that is pro¬ 
duced by the Boltzmann potential response to the paral¬ 
lel density gradients in the filamen t^^i^^’^^ . This can be 
achieved by isolating the component of 4>- according to 
Equation (1161) . The potential difference formed at high 
resistivities is thus demonstrated by the right hand plot 
of Figure [51 which plots the difference of 4>- between the 
mid-plane and sheath, from a i5x = 100, v\\ = 1000iz||q 
filament at the time of its peak radial velocity. For com¬ 
parison, the equivalent potential difference produced by 
a filament of the same (5x using the reference case resis¬ 
tivity, is illustrated in the left hand plot and is clearly 
negligible. 

Despite different physical mechanisms being dominant 
in determining the radial velocity of the filaments in 
the sheath current and resistive sheath current regimes 
(namely sheath resistivity and plasma resistivity respec¬ 
tively), the radial velocities produced in both regimes 
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FIG. 8. Potential difference formed between the mid-plane 
and sheath, ~ ~ ^00 filaments, 

at the time of their peak radial velocity. The left and right 
hand plot respectively the reference and enhanced 

resistivity v\\div = lOOOi^no cases. Each quantity is plotted at 
the time at which the filament’s characteristic radial velocity, 
Vf occurred. 



are proportional to {Tsheath + Tn) <5 _l The linear de¬ 
pendence on (Tg/ieat?! + Tji) is demonstrated in Figure 
[9] for (5 _l = 100 filaments. It is noted that a value of 
^sheath = 1/0.85 was used to plot this data. The high¬ 
est (Tsheath + Ty) data point deviates from this scaling 
because polarization currents were not negligible for this 
filament, and thus it is not strictly in the resistive sheath 
current regime. The transition between the sheath cur¬ 
rent and resistive sheath current regime therefore occurs 
at the point at which the Ty « T sheath, as this reflects 
the approximate point at which the sheath and plasma 
resistivities play an equal role in determining the fila¬ 
ment’s velocity. This is reflected in Figure [6l in that the 
filaments’ radial velocities only noticeably deviate from 
their reference case values once Fy « Tsheath ~ 1- This 
transition point and the observed velocity scalings in the 
two regimes are consistent with the predictions of Refer¬ 
ence [l3. 

Returning to the transition from the inertial regime to 
the sheath current or resistive sheath current regime, d* 
was measured quantitatively from the simulation data, 
by defining it to be the dj_ at which the maximum value 
oi Vf occurred, with a cubic spline interpolation used to 
determine its value as accurately as possible. The mea¬ 
sured locations of d* for each v\idiy simulated are plotted 
as blue crosses on the (5 _l - Fy diagram in Figure [in] to 
show the observed location of the boundary of the inertial 
regime. For reference, the horizontal dashed lines mark 
the values of F y corresponding to each of the values of 
i^Wdiv simulated, with the colors of the lines matching the 
colors of the markers used for their associated dataset in 
Figure m Also plotted using a dotted line is the location 
of the transition between the sheath current and resis¬ 
tive sheath current regimes. Moreover, the analytical es- 
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FIG. 9. Dependence of a = 100 filament’s characteristic 
radial velocity, Vf, on T sheath + Fy . With the exception of the 
highest Fsjteiith+ r|| point for each data series, in which polar¬ 
ization currents were not negligible, a linear dependence of u/ 
on r sheath +ry is displayed. Fy was scaled by a factor of 0.75 
when plotting the uniform izy data series, to account for the 
fact that the effective resistance to parallel currents traveling 
to the sheath from the mid-plane at 2 = 0 is approximately 
double that from the front of the hlament at 2 « F||/2. 
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timates for <5* from References]^ an dM which are stated 
in Equation dm, are plotted using solid black lines. For 
Fy < 1, the simulations’ d* remains constant around i5*o 
and is insensitive to F y and thus good agreement is found 
with the analytical predictions. For IT > 1, qualita¬ 
tive agreement is found with Reference I ITts prediction in 
that (5* increases as Fy rises. More quantitatively how¬ 
ever, the observed power law dependence in this region, 
^ Fy^^^'® Ri Ty®'^®, obtained from the two highest Fy 
data points and plotted using a blue dash-dot line, can 
be seen to have a weaker sc aling than the d* ~ 
scaling predicted in Reference [Tj. 


C. Increased Uniform Resistivity 

The effect of increasing the resistivity of the plasma 
uniformly in the SOL, rather than just in the region near¬ 
est to the targets, was also investigated to compare and 
contrast against the results in the preceding subsection. 
To enable a direct comparison, was increased by fac¬ 
tors such that the values of F y were approximately equal 
to those used in the sheath localized resistivity simula¬ 
tions. 

In terms of the effect onvf, for approximately the same 
value of Fy, increasing the resistivity uniformly through- 
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FIG. 10. Diagram of each filament regime location in Fn — 
space. The measured locations of the transition from the 
inertial regime, 5«, for different values of F|| are plotted using 
markers, whilst the analytical prediction from Reference [13 
for the location of 5* is plotted using a solid black line. The 
horizontal dashed lines mark the values of F y corresponding 
to those used in the sheath localized and uniform resistivity 
scans, with the colors of the lines matching the colors of the 
markers used for their associated dataset in Figures!^ or I 111 


out the domain produced very similar results to increas¬ 
ing it only in the last 25% nearest the target, as can be 
observed by comparing Figure [TT] to Figure El 

However, the actual balance of currents produced 
throughout the filament is subtly different. This can be 
seen by cross comparison between Figure [T^ which plots 
the divergence of current densities divided by n from the 
(5x = 28, = 2500jz||o, Fh =95 hlament, to Figure 

[71 Whilst in both cases the polarization current path is 
dominant in closing the diamagnetic current drive, the 
parallel currents can be seen to play less of a role in the 
uniform resistivity compared to the target localized re¬ 
sistivity case. This has two key effects. The first is that 
for comparable F ||, slightly higher radial velocities are at¬ 
tained in the inertial regime. For example, in Figure ITTl 
the zz|| = 2500 jz||o data series consistently attains higher 
radial velocities than v\\div = 10000zz||o data series in Fig¬ 
ure El for dx ^ 30. 

The second effect is that a differential radial velocity 
along the parallel direction is produced, in that the fila¬ 
ment moves faster in the radial direction at 2 : = 0 than 
it does at z = L||. This behavior is demonstrated by the 
left hand plot of Figure [T^l which plots n/ in a a: — z 
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FIG. 11. Dependence of the characteristic radial velocity, Vf, 
on its initial perpendicular length scale 5x for each of the 
values of izy used in the uniformly enhanced resistivity scan. 
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FIG. 12. Divergence of each current density divided by 
density, from the uniformly enhanced resistivity case 5x = 
28, I'll = 2500iZ||o filament at various perpendicular planes 
along the field line. The quantities are plotted at the time 
at which the filament’s characteristic radial velocity, Vf oc¬ 
curred. 
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FIG. 13. Comparison of the structure of the density pertur¬ 
bation, Uf , in an X — z plane through the middle of filaments 
using enhanced target localized resistivity u\\div = lOOOOi^no 
and uniformly enhanced resistivity = 2500r'||o. Both fil¬ 
aments were initialized with 5± = 12 and are shown at the 
time at which their characteristic radial velocities occurred. 
The values of resistivity are such that Fji is approximately 
equal in both cases. 


plane through the middle of the (5j_ = 12, ^|| = 2500iz||o 
filament at the time at which its characteristic velocity, 
Vf, occurs. For comparison the same quantity from the 
target localized resistivity v\\div = 10000iz||q, (5j_ = 12 
simulation is plotted in the right hand plot of the same 
figure. By introducing the drift plane radial velocity. 


= -. ( 19 ) 


oo oo 

J J Uf dxdy 


— OO —OO 


the extent to which the radial velocity varies along the 
field line can be quantitatively assessed. Figure [Ml plots 
this quantity against time at various positions along the 
field line in the region of density perturbation for the 
i5_l = 12 simulations shown in Figure [13] in addition to 
the reference case (5j_ = 12 simulation. It can be seen 
that in the uniformly enhanced resistivity case, the radial 
velocity at z = F||/2 is approximately a third slower than 
at z = 0 for most of the simulation. 

The differential radial velocity along the field line oc¬ 
curs in the enhanced uniform resistivity case because the 





FIG. 14. Drift plane radial velocity, Ua,, plotted against time 
at various positions along the field line, for 5± = 12 filaments. 
The top, middle and bottom plots respectively correspond to 
the reference resistivity case, the uniformly enhanced resis¬ 
tivity zz|| = 2500iZ||o case and the enhanced sheath localized 
resistivity Vf^div = 10000iZ||o case. A differential radial velocity 
along the field line is exhibited in the izy = 2500 :z||q case. 


parallel currents, which are suppressed throughout the 
domain, are not able to balance the polarization currents 
in the region T||/2 < z < 3L||/4, where the diamag¬ 
netic currents are reduced or negligible. The polariza¬ 
tion currents therefore can only develop to match the lo¬ 
cal diamagnetic current drive, meaning that they are not 
constant along z. Consequently (j) and Vx develop larger 
values at z = 0 than at z = T||/2, where the density gra¬ 
dients are smaller. If the resistivity is high enough, (j) is 
determined locally on each drift plane and the dynamics 
of the hlament are effectively decoupled along the field 
line. In contrast, in the enhanced target localized resis¬ 
tivity case, the parallel currents are able to fulfill the role 
of balancing the polarization currents in the region where 
the diamagnetic currents are absent, and so Jpoi, 4> and 
Vx are approximately constant from z = 0 to z = 3L||/4. 

Regarding the behavior of under uniformly in¬ 
creased resistivity, a very similar trend was displayed to 
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what was found using an enhanced target localized re¬ 
sistivity. This can be observed in Figure (TUI which plots 
using red circles the measured locations of from this 
uniform resistivity series of simulations. The measured 
power law dependence, ^ Tji ' ' ~ Fn ' in the re¬ 
gion Fji > 1 is approximately the same as in the sheath 
localized resistivity case and is again weaker than Refer¬ 
ence [i3’s prediction. 

V. CONCLUSIONS 

In this work the influence of enhanced parallel resistiv¬ 
ity on the dynamics of SOL filaments has been studied 
using 3D simulations. Motivated by the expectation of 
lower temperatures, high neutral densities and possibility 
of detachment in the divertor region, the resistivity was 
increased only in the last quarter of the domain nearest 
the targets. Increasing the resistivity lead to a suppres¬ 
sion of parallel currents, a corresponding enhancement 
of polarization currents, and the development of a po¬ 
tential difference along the field line. These intrinsically 
3D effects meant that filaments attained higher radial ve¬ 
locities at enhanced resistivities. In particular, filaments 
with a large perpendicular length scale, (5_l, experienced 
the greatest increase in radial velocity, because at low 
resistivity these filaments were sheath current limited, 
meaning that their diamagnetic currents were predomi¬ 
nantly closed via parallel currents. In contrast, polariza¬ 
tion currents were dominant in ensuring current continu¬ 
ity for the smallest (5_l filaments at low resistivities and 
so these filaments only experienced a modest increase in 
their radial velocities at higher resistivities. 

More specifically, one mechanism by which greater ra¬ 
dial velocities were produced at higher resistivities, is 
that filaments that were in the sheath current regime 
transitioned into the inertial regime and so velocities 
scaled like 5 _l up until larger values of (5 _l . The critical 
(5 _l at which this transition occurs, (5*, thus increased with 
resistivity and its dependence on the total resistance to 
parallel currents between the mid-plane and sheath en¬ 
trance, Fji, was measured to be approximately <5* oc Ljj-^, 
which is marginally weaker than that predicted by Ref¬ 
erence [B Enhanced radial velocities were also observed 
in filaments that were sufficiently large (5 _l ^ 6^,), such 
that the parallel current path remained dominant over 
the polarization current path even at enhanced resistiv¬ 
ities. The mechanism for these resistive sheath current 
regime filaments was that the resistance of the plasma 
was sufficient to introduce a potential difference between 
downstream at the sheath entrance and further upstream 
in the region of the filament density perturbation, such 
that for the same amount of current to flow into the 
sheath, larger potentials were formed upstream at higher 
resistivity, corresponding to faster radial velocities. 

Investigations were also carried out in which the resis¬ 
tivity was increased uniformly throughout the domain. 
For the same value of Fji, marginally faster radial veloci¬ 


ties were produced in the uniform resistivity case, as the 
parallel currents were more effectively suppressed. The 
biggest difference with respect to the target localized re¬ 
sistivity simulations however was that filaments exhib¬ 
ited a differential radial velocity along the field line in 
the uniformly enhanced resistivity case, moving radially 
faster at the mid-plane than further downstream. This 
demonstrated that enhanced resistivity can decouple the 
dynamics of filaments along the parallel direction. 

A limitation of the simulations presented in this paper 
is that since the model used assumes isothermal electrons 
and neglects neutral physics, the target localized resistiv¬ 
ity was arbitrarily increased rather than self-consistently 
calculated. Moreover the simulations have neglected the 
influence of enhanced magnetic shear around the X-point 
region, which could provide alternative current paths for 
the diamagnetic currents to be closed. Furthermore, the 
assumption of cold ions is poorly justified in the SOL, 
where typically Ti > T^. Previous 2D gyrofluid model 
simulations2ii2^ have demonstrated that the inclusion of 
finite Larmor radius effects increase the coherence of fila¬ 
ments as they move radially outwards and induces the fil¬ 
ament to move in the poloidal direction. Therefore the in¬ 
clusion of electron temperature dynamics, hot ions, neu¬ 
tral physics and magnetic geometry effects to this model 
would be useful additions for future research. Never¬ 
theless, this work has demonstrated the mechanisms by 
which enhanced divertor resistivities may produce faster 
radial filament velocities. 
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Appendix: Estimate of the magnitude of u^di-v at low 
divertor temperatures 

This appendix provides an order of magnitude esti¬ 
mate of the electron temperature in the divertor region, 
that may be required to produce the values of vwdiv 
used in the target localized resistivity study in Section 
IIVBI The definition of given in Section m is based 

upon electron-ion collisions, and is such that ocTe . 
However, at very low temperatures, collisions between 
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electrons and neutrals may become important in the di¬ 
vertor region, and so more generally, can be definedi^ 
as: 

'^ll “ ^if + ■ (^-1) 

Here = i^eio/(l-96Hi) is the normalized electron-ion 
collisionality (given as the definition of in Section HIl) . 
and = Vg^n! (l-960ine/no) is the normalized electron- 
neutral collisionality where Ug is the density of electrons. 
Furthermore Vgn = rin (cu) where n„ is the density of 
neutral atoms, v is the velocity of electrons, a is the 
cross section for collisions between electrons and neu¬ 
trals (and in principle is a function of v), whilst (■) 
denotes averaging over all velocities in the (assumed) 
Maxwellian distribution function. The densities of neu¬ 
tral deuterium atoms and electrons (or deuterium ions) 
at a given temperature were estimated using the Saha 
equation^, alongside the assumption that Ug + Un = no. 
It is noted that the Saha equation assumes the plasma 
and neutral gas to be in thermal equilibrium, which may 
not be a valid assumption for edge plasmas. Moreover, 
its use implies an equilibrium between ionization and re¬ 
combination processes, which may not occur because re¬ 
combination is a relatively slow process compared to the 
typical timescales of fluctuations in the divertor. Values 
of a for elastic collisions between electrons and hydrogen 
atoms, obtained from Referencels^, were used. These cal¬ 
culations arguably provide a conservative estimate of the 
resistivity in the divertor region, as collisions with neu¬ 
tral particles other than deuterium, that may be present 
due to sputtering or impurity seeding, have not been in¬ 
cluded. Furthermore, anomalous resistivity effects have 
also been neglected. 

Figure dHl thus shows the estimated relative increase of 
the normalized collisionality in the divertor region, 

(and each of its constituent terms, and J^y^^^) over 

i^llo as the temperature in the divertor, is decreased. 
Moreover, for reference, the corresponding values of 
estimated to produce the values of including and 
excluding neutral collisions are given in Table HIl These 
calculations indicate that at temperatures around O.SeV, 
electron collisions with neutrals may dominate over elec¬ 
tron collisions with ions. 

It is emphasized however, that the calculations pre¬ 
sented in this Appendix are based upon a number of as¬ 
sumptions that may not be well justihed in the SOL near 
the sheath, and so should only be used to give an order of 
magnitude indication of how the resistivity in the diver¬ 
tor region may depend on Tg. Moreover, it is important 
to note that J^y" dominates at low temperatures because 
the ratio n„/ne becomes very large (> 1000) and so the 
plasma is estimated to be very weakly ionized. It is un¬ 
clear whether such a weakly ionized plasma is achieved 
experimentally in the divertor region when detachment 
occurs. 
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FIG. 15. Estimate of relative magnitude of the normalized 
collisionality in the divertor region, (and each of its con¬ 
stituent terms, and compared to izyQ as a function 

of temperature in the divertor, Tf™. 

TABLE II. Estimated divertor temperature required to pro¬ 
duce the values of used in Section flV Bl 


^\\div 

Til 

Tf " required 
excluding neutral 
effects (eV) 

rf " required 
including neutral 
effects (eV) 

V\\0 

0.038 

40 

40 

10!Z||o 

0.12 

8.6 

8.6 

O 

O 

o 

0.98 

1.9 

1.9 

1000p||o 

9.6 

0.40 

0.45 

10000p||o 

95 

0.086 

0.37 

100000i/||o 

950 

0.019 

0.32 
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